#!usr/bin/perl -w

use strict;
use FindBin qw($Bin $Script);
die "fa.list format\nxx1\txx1.fa\tcorr-name1\nxx2\txx2.fa\tcorr-name2\n..\nusage:perl $0 <fa.list> <muscle.list> \n" unless @ARGV == 2;

open A,$ARGV[0];    open B,$ARGV[1];
my $num=0;my %h; my %gene; my %gene_id; my @sam;
while(<A>){
	chomp;my @a=split;  $num++;  	open IN,$a[1];	
	if($a[2] ne "WHLR"){
		push @sam,$a[2];
	}
	$/=">";    <IN>;	
	while(my $inf=<IN>){
		chomp $inf;my @b=split /\n+/,$inf;my @c=split /\s+/,$b[0];
		if(@a == 3){			$gene_id{$c[0]}=$a[2];		}
		else{					$gene_id{$c[0]}=$a[0];		}
	}
	close IN;
	$/="\n";
}
close A;	
while(<B>){
	chomp;my @a=split;
	open IN,$a[1];$/=">";<IN>;my %part;
	my $whl;
	while(my $inf=<IN>){
		chomp $inf;my @b=split /\n+/,$inf;	my $seq=join "",@b[1..$#b];	my @c=split /\s+/,$b[0];
		if($gene_id{$c[0]} eq "WHLR"){
			$whl=$seq;
		}
		else{
			$part{$gene_id{$c[0]}}=$seq;
		}
	}
	$/="\n";
	foreach my $k1(@sam){
		my $all=0;	my $sam=0;
		my @sub1=split //,$part{$k1};
		my @sub2=split //,$whl;
		#print "@sub1\n@sub2\n";
		foreach my $k2(0..$#sub1){
			next if($sub1[$k2] eq "-" || $sub2[$k2] eq "-");
			$all++;
			if($sub1[$k2] eq $sub2[$k2]){
				$sam++;
			}
		}
		if($all == 0){
			print STDERR "$a[0]\t$k1\n";
			$h{$a[0]}{$k1}=0.5;
		}
		else{
			$h{$a[0]}{$k1}=$sam/$all;
		}
	}
}
close B;
my $sam=join "\t",@sam;
print "ID\t$sam\n";
foreach my $k1(sort keys %h){
	print "$k1";
	foreach my $k2(@sam){
		print "\t$h{$k1}{$k2}";
	}
	print "\n";
}

